%%matlabpool close
%matlabpool open 

%% Computational Vision for EM


clear all;
close all;
clc;
folder = 'Urban 12'; %% Define which folder to test
file_List = ls(folder);
num_Pic = size(ls(folder));
num_Pic = num_Pic(1) - 3;

fprintf('Folder Set to "%s"\n',folder);


%% Image sub set Classify
MaxT_BetweenFrame = 3600; %define the maximum time between two frame 

nb_vector=1;
Set_Array = 0;


ImF=imfinfo(strcat(folder,'/',file_List(3,:)));

h=str2num(ImF.DateTime(12:13));
m=str2num(ImF.DateTime(15:16));
s=str2num(ImF.DateTime(18:19));
T=s+60*m+60*60*h;


Set_Array(nb_vector,1) = 3;
Set_Array(nb_vector,2) = 0;
if strcmp(ImF.DigitalCamera.SceneCaptureType,'Standard')
    Set_Array(nb_vector,3) = 'd';
else
    Set_Array(nb_vector,3) = 'n';
end
    
for l=2:num_Pic

    ImF=imfinfo(strcat(folder,'/',file_List(l+2,:)));
    
    h=str2num(ImF.DateTime(12:13));
    m=str2num(ImF.DateTime(15:16));
    s=str2num(ImF.DateTime(18:19));
    T2=s+60*m+60*60*h;
    
    %if abs(T2-T)>MaxT_BetweenFrame
    
    if strcmp(ImF.DigitalCamera.SceneCaptureType,'Standard')
        this_set = 'd';
    else
        this_set = 'n';
    end
    
    if this_set ~= Set_Array(nb_vector,3) || abs(T2-T)>MaxT_BetweenFrame
        nb_vector=nb_vector+1;
        Set_Array(nb_vector,1) = l+2;
        Set_Array(nb_vector,2) = 0;
        if strcmp(ImF.DigitalCamera.SceneCaptureType,'Standard')
            Set_Array(nb_vector,3) = 'd';
        else
            Set_Array(nb_vector,3) = 'n';
        end
    else
        Set_Array(nb_vector,2)=Set_Array(nb_vector,2)+1;
    end
    T=T2;
end

clear T T2 l h m s ImF this_set;

%%
for g = 1:nb_vector
    Load_count =  Set_Array(g,2)+1;
    Start_count =  Set_Array(g,1)-1;
    
    if Load_count== 0
        continue;
    end

%% Image Preprocessing

fprintf('Loading Images and Calculate Edges.....');

clear im_back_ground im_back_edge;
tmp_im = rgb2gray(imread(strcat(folder,'/',file_List(3,:))));
pic_size = size(tmp_im);
if pic_size(2) == 1024
    tmp_im = tmp_im(17:752,:); %turkey
elseif pic_size(2) == 1920
    tmp_im = tmp_im(33:1048,:);
end

tmp_im = tmp_im(1:3:end,1:3:end,:);
tmp_size = size(tmp_im);
im_back_ground = zeros(tmp_size(1),tmp_size(2),Load_count);
im_back_edge = zeros(tmp_size(1),tmp_size(2),Load_count);
for i = 1:Load_count
    tmp_im = rgb2gray(imread(strcat(folder,'/',file_List(Start_count+i,:))));
    if pic_size(2) == 1024
        tmp_im = tmp_im(17:752,:); %turkey
    elseif pic_size(2) == 1920
        tmp_im = tmp_im(33:1048,:);
    end
	im_back_ground(:,:,i) = tmp_im(1:3:end,1:3:end,:);
	im_back_edge(:,:,i) = edge(im_back_ground(:,:,i),'canny');
end

fprintf('Finished\n');

fprintf('Estimating Image Edge.....');

edge_sum = sum(im_back_edge,3);
est_e_mean = sum(sum(edge_sum))/sum(sum(edge_sum~=0));
est_edge = edge_sum > est_e_mean;

fprintf('Finished\n');

fprintf('Calculating Images Mean and Var.....');

clear tmp_im tmp_size;

est_para = Background_Estimate('GAUSSIAN',im_back_ground);

fprintf('Finished\n');


fprintf('Creating Hollow Image.....');

warning off all
clear new_background;
need_b = 0;
im_size = size(im_back_ground(:,:,1));
new_background = est_para(:,:,1);
glb_var = mean(reshape(sqrt(est_para(:,:,2)),1,[]));
n_classes = zeros(im_size(1),im_size(2),2);
for i = 2:im_size(1)-1
    for j = 2:im_size(2)-1
        if est_para(i,j,2) > (glb_var^2)
            new_background(i,j) = 0;
            tmp = im_back_ground(i,j,:);
            tmp = tmp(:);
            [IDX, C] = kmeans(tmp, 2,'Start','uniform','EmptyAction','drop');
            n_classes(i,j,:) = C;
            %[N,X]  = hist(tmp,255);
            %n_classes(i,j,:) = separate_class(N,glb_var);
        end
        if new_background(i,j) == 0
            need_b = need_b +1;
        end
    end
end
warning on all
clear tmp tmp_class tmp_xout im_back_edge;
bak_new_background = new_background;

fprintf('Finished\n');


%%
warning off all
new_background = bak_new_background;
need_b_ct = need_b;
prev_b_ct = need_b_ct;

fprintf('Estimate Background Image.....');
p = 0;

m_min_c = 4;
cur_var = 1;
tmp_cur_var = 1;
rev_pro = 0;
Process_back = 1;
i_ = (cur_var+1):im_size(1)-(cur_var);
j_ = (cur_var+1):im_size(2)-(cur_var);
f_f = 0;
FillUp = 0;

%imwrite(uint8(new_background),strcat('.\BackGround est\DEMO\Night\',int2str(need_b_ct),'.JPG'),'jpg');

while(Process_back == 1 && need_b_ct >= 0)
    Process_back = 0;
    if f_f == 0
        f_f = 1;
        i_ = fliplr(i_);
    else
        f_f = 0;
        j_ = fliplr(j_);
    end
    for i = i_
        for j = j_
            if(new_background(i,j) == 0) 
                clear tmp_background;
                if(est_edge(i,j)== 0)
                    MeEdge = 0;
                else
                    MeEdge = 1;
                end
                if rev_pro == 1 && MeEdge == 1
                    continue;
                end
                tmp_background = zeros(1,1);
                m_count = 1;
                for i_b = [i-cur_var:i+cur_var]
                    for j_b = [j-cur_var:j+cur_var]
                        if (i == i_b && j == j_b) 
                            continue;
                        end
                        if( new_background(i_b,j_b) ~= 0 )
                            b_diff = sqrt(glb_var);
                            for n_cls = 1:2
                                if(n_classes(i,j,n_cls)~= 0 && (est_edge(i_b,j_b) == MeEdge || FillUp))
                                    tmp_diff = abs(int32(n_classes(i,j,n_cls))- int32(new_background(i_b,j_b)));
                                    if tmp_diff < b_diff
                                    	tmp_background(m_count) = n_classes(i,j,n_cls); 
                                    end
                                end
                            end
                            if(tmp_background(end) ~= 0)
                                m_count = m_count + 1;
                                tmp_background(m_count) = 0;
                            end
                        end
                    end
                end
                if (m_count > m_min_c)
                    if mode(tmp_background) ~= 0
                        need_b_ct = need_b_ct - 1
                        Process_back = 1;
                        new_background(i,j) = mode(tmp_background);
                    end
                end
            end
        end
    end
    
        if ( (prev_b_ct - need_b_ct)/prev_b_ct < 0.3 )
            if need_b_ct > 0
                Lky_size = 3;
                Lky_the = (Lky_size*2 + 1)^2;
                [Lky_IDX, Lky_C] = kmeans(nonzeros(new_background),2,'Start','uniform');
                if Set_Array(g,3) == 'n'
                    Lky_int = 0;
                else
                    Lky_int = median(nonzeros(new_background));
                end
                for Lky_i = (Lky_size+1):im_size(1)-(Lky_size)
                    for Lky_j = (Lky_size+1):im_size(2)-(Lky_size)
                        Lky_i_b = Lky_i-Lky_size:Lky_i+Lky_size;
                        Lky_j_b = Lky_j-Lky_size:Lky_j+Lky_size;
                        if new_background(Lky_i,Lky_j)~=0 || sum(sum(new_background(Lky_i_b,Lky_j_b)==0)) > Lky_the*0.9
                            continue;
                        end
                        Lky_tmp_background = zeros(1,1);
                        Lky_m_count = 1;
                        for Lky_i_b = Lky_i-Lky_size:Lky_i+Lky_size;
                            for Lky_j_b = Lky_j-Lky_size:Lky_j+Lky_size;
                                Lky_b_diff = abs((Lky_C(2)-Lky_C(1))/2);
                                if FillUp == 1
                                    Lky_b_diff = glb_var * tmp_cur_var;
                                end
                                for Lky_n_cls = 1:2
                                    if(n_classes(Lky_i,Lky_j,Lky_n_cls) ~= 0 )
                                        Lky_tmp_diff = abs(int32(n_classes(Lky_i,Lky_j,Lky_n_cls))- int32(Lky_int));
                                        if Lky_tmp_diff < Lky_b_diff
                                            Lky_tmp_background(Lky_m_count) = n_classes(Lky_i,Lky_j,Lky_n_cls);
                                        end
                                    end
                                end
                                if(Lky_tmp_background(end) ~= 0)
                                    Lky_m_count = Lky_m_count + 1;
                                    Lky_tmp_background(Lky_m_count) = 0;
                                end
                            end
                        end
                        if (Lky_m_count > Lky_the*0.9)
                            if mode(Lky_tmp_background) ~= 0
                                new_background(Lky_i,Lky_j) = mode(Lky_tmp_background);
                                need_b_ct = need_b_ct - 1
                                Process_back = 1;
                            end
                        end
                    end
                end
                clear Lky_i Lky_j Lky_size Lky_i_b Lky_j_b Lky_int Lky_tmp_background Lky_m_count Lky_n_cls Lky_IDX Lky_C;
            end
        end
    if((prev_b_ct - need_b_ct)/prev_b_ct < 0.1 )
        tmp_cur_var = tmp_cur_var + 1;
        if glb_var*tmp_cur_var < 255
            Process_back = 1;
        end
    end
    if(Process_back == 0 || (prev_b_ct - need_b_ct)/prev_b_ct < 0.01)
        if FillUp ~= 1
            FillUp = 1;
            Process_back = 1;
        end
    end
    prev_b_ct = need_b_ct;
end

figure
imshow(new_background,[]);
fprintf('Created\n');



%%

clear im_gray;
parfor i = 1:Load_count
    im_gray(:,:,i) = Background_Subtruction('BASIC',im_back_ground(:,:,i),new_background,0);
end

for k = 1:Load_count
    imwrite(im_gray(:,:,k),strcat('.\TMP\',int2str(Start_count+k),'.JPG'),'jpg');
end

%% This Part is for removing the blocks that appeared in the image too
%% often ( ususlly cause by error background estimation)

%{
fprintf('Image Subtruaction Complete\n');

im_gray_bak = im_gray;
im_gray = im_gray_bak;

fprintf('Estimating Foreground.....');
est_Gray = Background_Estimate('GAUSSIAN',im_gray);
fprintf('Finished\n');

est_size = size(est_Gray);
est_V_mean = sum(sum(est_Gray(:,:,2)))/sum(sum(est_Gray(:,:,2)~=0));
tmp_V_v = est_Gray(:,:,2) - est_V_mean;
tmp_V_v = tmp_V_v.^2;
est_V_v = sum(sum(tmp_V_v))/sum(sum(tmp_V_v~=0));
clear tmp_V_v;

fprintf('Removing Possible Error.....');

for i = 1:est_size(1)
    for j = 1:est_size(2)
        if est_Gray(i,j,2) < est_V_mean + sqrt(est_V_v)
           im_gray(i,j,:) = im_gray(i,j,:) - mean(im_gray(i,j,:));
           for k = 1:Load_count
               if im_gray(i,j,k) < 0
                   im_gray(i,j,k) = 0;
               end
           end
        end
    end
end

for k = 1:Load_count
    imwrite(im_gray(:,:,k),strcat('.\BackGround est\Est Set\tmp\B Sub with ForeGround Est\',int2str(k),'.JPG'),'jpg');
end

fprintf('Finished\n');

fprintf('Threadshold the Images.....\n');
im_gray(:,:,:) = abs(im_gray(:,:,:));
for k = 1:Load_count
    est_g_mean = sum(sum(im_gray(:,:,k)))/sum(sum(im_gray(:,:,k)~=0));
    tmp_g_var = im_gray(:,:,k) - est_g_mean;
    tmp_g_var = tmp_g_var.^2;
    est_g_var = sum(sum(tmp_g_var))/sum(sum(tmp_g_var~=0));
    im_resut = im_gray > est_g_mean + sqrt(est_g_var)*.5;
end
fprintf('Finished\n');

fprintf('Result will be shown\n');

for k = 1:Load_count
    imwrite(im_resut(:,:,k),strcat('.\BackGround est\Est Set\tmp\B Sub with ForeGround Est Thre\',int2str(k),'.JPG'),'jpg');
end


for k = 1:Load_count
    figure
    imshow(im_gray(:,:,k),[]);
end
%}

end


